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Abstract 

An asymptotic treatment of thin accretion disks, introduced by Kluzniak & Kita (2000) for a steady-state disk flow, 
is extended to a time-dependent problem. Transient growth of axisymmetric disturbances is analytically shown to 
occur on the global disk scale. The implications of this result on the theory of hydrodynamical thin accretion disks, 
as well as future prospects, are discussed. 
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1. Introduction 

Thin accretion disks (AD) form whenever a suffi- 
ciently cool gas, endowed with a significant amount 
of angular momentum, is gravitationally attracted 
towards a relatively compact object. This situation 
is quite common in astrophysics and therefore the 
observational and theoretical study of accretion 
disks has been quite intensive. The necessary condi- 
tion for accretion to take place is that angular mo- 
mentum be extracted from the fluid swirling around 
the central object. To achieve accretion rates that 
are consistent with observations, the physical mech- 
anism for such angular momentum transport must 
be more efficient (by many orders of magnitude) 
than just the one resulting from torques caused by 
microscopic viscosity. 

Already at the outset, when AD were theoret- 
ically proposed (Prendergast & Burbidge 1968, 
Pringle & Rees 1972), the enhanced effective vis- 
cosity was postulated to result from turbulence 
and was parametrized, using mixing length theory 
or a similar scheme, as detailed theoretical under- 
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standing of turbulence and the transition to it was 
lacking then (a situation that is still with us nowa- 
days). The parametrization of the effective viscosity 
in disks in terms of a single parameter - a, intro- 
duced by Shakura & Sunyaev (1973), proved itself 
to be the most fruitful, giving rise to successful 
interpretations of many observational results (see 
Lin & Papaloizou 1996, Frank, King & Raine 2002, 
for reviews). In some cases (e.g. dwarf nova mod- 
els), however, more complex (and cumbersome) 
prescriptions for the viscosity had to be employed. 
However, until the early 1990's the question what is 
the physical origin of turbulence (or, more precisely, 
the anomalous angular momentum transport) in 
AD has essentially remained unanswered. No def- 
inite linear instability has been identified for thin 
disks in which the swirling flow consist of Keplerian 
shear. The magneto- rotational instability (MRI), 
originally found (Velikhov 1959, Chandrasekhar 
1960) for magnetic Taylor-Couette, (i.e., cylindri- 
cal) flows, has been shown by Balbus & Hawley 
(1991) to also operate in cylindrical AD, when the 
fluid is electrically conducting and for not too large 
initial magnetic fields. 

It has thus become the operating paradigm in 
the astrophysical community that purely hydrody- 
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namical turbulence in Keplerian disks is altogether 
ruled out and since the problem of angular momen- 
tum transport in these objects relies on MRI driven 
magneto-hydrodynamical (MHD) turbulence, ex- 
tensive numerical calculations are necessarily the 
main tool of research on this problem (sec Balbus 
& Hawley 1998 and Balbus 2003, for reviews and 
references). Nonetheless, efforts to find a purely 
hydrodynamical transition to turbulent activity 
in thin AD have continued until the present time, 
more than a decade after Balbus , Hawley & Stone 
(1996) appeared to settle the matter. These have 
largely been motivated by the fact that a purely hy- 
drodynamical AD flow has a 'microscopic' Reynolds 
number (Re) of the order of 10 14 or so, a quite 
amazing setting for a laminar flow. 

Among the ideas that have been put forward in 
this context some are based on the application to 
thin disk flows of the viewpoint, familiar to the fluid- 
dynamical community, that transient dynamics in- 
duced by perturbations, i.e. strong transient growth 
(TG) in linearly stable shear flows may play an im- 
portant role in nonlinearly shaping the final dynam- 
ical state. TG is possible because the relevant lin- 
ear operator governing the behavior of infinitesi- 
mal perturbations in these flows is non-normal and 
thus the usual modal approach essentially fails (see, 
e.g., Grossman 2000, Schmid & Hcnnigson 2001, and 
Criminale, Jackson & Joslin 2003) 

TG of disturbances has been discussed in the as- 
trophysical literature within two quite different (but 
related) contexts: 

• Local disturbances experiencing large enough 
TG that can possibly trigger a subcritical nonlin- 
ear transition into turbulence in a linearly stable 
shear flow, i.e. via the so-called bypass transi- 
tion (e.g., Chagelishvili et al. 2003, Tevzadze 
et al. 2003, Yecko 2004, Afshordi et al. 2005, 
Mukhopadhyay et al. 2005) 

• Perturbations experiencing significant TG that 
can be excited either by some external agent, per- 
haps as secondary flows on a pre-existing 3D tur- 
bulence itself, however weak, may give rise to in- 
tense global dynamical activity (e.g. Ioannaou & 
Kakouris 2001, and see below). 

Other scenarios for inducing hydrodynamical ac- 
tivity in AD (some of them directly or indirectly 
related to TG as well) have also been proposed. 
Among them are those invoking baroclinic insta- 
bilities (e.g., Klahr & Bodcnhcimcr 2003), strato- 
rotational instabilities (e.g., Dubrullc et al. 2005, 
Umurhan 2006), formation and long sustenance of 



vortices and/or waves (e.g., Godon & Livio 1999, 
Bracco et al. 1999, Li et al. 2001, Umurhan & Regev 
2004, , Barranco & Marcus 2005, Petersen, Stew- 
art & Julien 2007, Bodo et al. 2007, Lithwick 2007). 
Even though no definite and undisputed hydrody- 
namical mechanism, that can give rise to sufficient 
angular momentum transport in AD, has been iden- 
tified so far, significant efforts along these lines are 
continuously being made by several groups of re- 
searchers. 

It appears that continued study of purely hydro- 
dynamical processes in disks remains viable and 
worthwhile, in particular in view of the difficulties 
with MRI driven angular momentum transport in 
AD, which have recently been pointed out. Ques- 
tions of insufficcnt numerical resolution in MHD 
disk simulations have been convincingly raised (Pc- 
sah, Chen & Psaltis 2007, Fromang & Papaloizou 
2007). The decrease of transport with decreasing 
magnetic Prandtl number (Pm) for various setups 
and boundary conditions and, in particular, the 
vanishing of this transport for Pm <C 1 (which is 
the case in many instances of astrophysical AD) 
has been demonstrated (Umurhan, Mcnou & Regev 
2007, Lesur & Longaretti 2007, Umurhan, Regev & 
Menou 2007, Fromang et al. 2007). Finally, serious 
doubts as to the viability of the local shearing box 
approximation (Goldrcich & Lynden-Bell 1965) to 
the numerical study of accretion disk MHD turbu- 
lence (Coppi & Keyes 2003, King, Pringle & Livio 
2007, Shu et al. 2007, Regev & Umurhan 2007) have 
been raised. 

The purpose of this contribution is to draw at- 
tention to the possibility that a thin AD, in which 
a sub-critical transition to hydrodynamical turbu- 
lence occurs (for a sufficiently high Re) large tran- 
sient amplification of global disturbances may in- 
duce recurrent or even persistent secondary flow 
activity. Lesur & Longaretti (2005) have explicitly 
demonstrated, using high resolution 3D numerical 
simulations, that such sub-critical transition to tur- 
bulence does appear in a model flow having Keple- 
rian shear characteristics, but at a very high Re (see 
also Rincon et al. 2006). However, it appears that 
the efficiency of turbulent transport in these flows is 
insufficient for astrophysical purposes (i.e. AD). A 
similar conclusion can perhaps be drawn from the 
recent experimental study of Ji et al. (2006). The 
very small value of the effective a (measuring the 
angular momentum transport) in these flows is, as 
we shall show, the key necessary ingredient for the 
excitation of vigorous global (transient) secondary 
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flows, atop the weak turbulent state. It is quite con- 
ceivable that the large (by orders of magnitude in 
the disturbance energy) TG may (nonlinearly) give 
rise to persistent dynamical activity, or at least be 
recurrently re-excited by external perturbations (see 
Ionnaou & Kakouris 2001, who advocated the lat- 
ter possibility). It however remains to be shown, of 
course, that in the ultimate state angular momen- 
tum transport is appropriate for AD, i.e. the effec- 
tive a acquires a high enough value. 

The primary tool utilized here in order to facili- 
tate an analytical treatment is the asymptotic ex- 
pansion, where the dependent variables and gov- 
erning equations are expanded in powers of a small 
quantity (here the measure of the disk's 'thickness', 
e - see below) . Exposing the resulting mathematical 
system (an initial value problem) to a set of initial 
conditions, in which the extreme geometry of the 
disk structure is taken into account (i.e. e <C 1), we 
aimed at obtaining an analytically treatable prob- 
lem. To achieve this goal we also assumed a poly- 
tropic relation between the pressure and density. 
The advantage of such an approach is obvious - the 
treatment can be essentially analytical and the re- 
sponsible physical effects leading to any interesting 
dynamics may be transparently traced. The study 
reported on here is limited to axi-symmetric distur- 
bances and the question of its ultimate development 
still remains open. We are now in the process of gen- 
eralizing it to fully 3-D disturbances and it is clear 
that the present analytic work should be ultimately 
complemented with detailed and uncompromising 
3-D numerical calculations with a proper treatment 
of energy generation and transfer. 

2. Asymptotic formulation of the dynamics 
of a thin polytropic accretion disk 

We wish to investigate the dynamical (i.e. time 
dependent) behavior of thin, axisymmetric AD, al- 
lowing for vertical structure. Viscous flow is invoked, 
with a standard viscosity prescription, that is we ac- 
tually assume that an effective viscosity producing 
mechanism is already operative. The mean flow is 
thus assumed to consist of Keplerian thin disk flow, 
which has already undergone a sub-critical tran- 
sition, like the one found by Lesur & Longaretti 
(2005), for example. As mentioned above, these au- 
thors found that the value of the effective a scales 
with the value of the sub-critical transition Reynolds 
number, Re 9 , asa~ 1/Re 9 and since they obtained 



that the value of Re g for Keplerian rotation flows 
is extremely high, they concluded that the a values 
that can be expected are much smaller than those 
necessary for angular momentum transport in as- 
trophysical AD. He we shall examine the dynamical 
behavior for different values of a and, as we shall 
see, the TG we find will be anti-correlated with the 
value of a of the mean steady flow. 

The key observation is that when on writes the 
equation of fluid dynamics (e.g., Tassoul 1978) and 
makes them non-dimensional by scaling the vari- 
ables with their typical values, a non-dimensional 
parameter e appears, multiplying (by its value in 
various powers) different terms. These (axisymmet- 
ric) equations in cylindrical coordinates (r, z, 4>) 
read: 



dtp + -d r {rpu) + d z (pv) = 0, 
r 



(1) 



ed t u + e 2 ud r u + evd z u — fl 2 r = — - ( 1 + e 2 — 



+ -d z {r)d z u) + — 
P P 



pd r vj + d z (rjd r v) — -d r (rjd z v) 



p y r z r 3 Lr 



d t v + ed r v + vd z v = I 1 + e 2 ^- ] - d z w 



+ 



(2) 



^d z {r]d z v) + ^- \d r (i]rd z u) - ^d z [rfd r (ru)]^ 



(3) 



+ e — d r (rird r v) . 
pr 



+-d z (f]d z n) +^d r ( V r 3 d r n) . (4) 

p pr° 

In these equations u, v and O are the horizontal, ver- 
tical and angular velocity respectively, w = nc 2 s = 
(n + X)p x l n (c s is the sound speed and the non- 
dimensional polytropic relation P = p 1+1 / n , with 
the general polytropic index n, has been used). Also, 
77 is the dynamic viscosity, which can be expressed 
using the standard Shakura-Sunyaev viscosity pre- 
scription 
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with the non-dimensional Keplerian angular veloc- 
ity substituted as £1k = r~ 3 / 2 . 

The non-dimensional parameter e = c s /(f2r) = 
h/r (the quantities denoted by tilde are the corre- 
sponding dimensional quantities, evaluated at a typ- 
ical position in the AD, with h being the disk height) 
is very small for thin disks (that is, those assumed 
to be able to cool efficiently, so that the rotational 
velocity is highly supersonic). 

An asymptotic approach of the kind used here 
was introduced for the first time to the study of thin 
viscous AD by Regev (1983), in the context of AD 
boundary layers and later was developed and used 
in a remarkable analytical work of Kluzniak & Kita 
(2000), (hereafter KK) who solved for the steady 
structure of a polytropic viscous axisymmetric disk. 
This study revealed the presence of a steady merid- 
ional flow pattern with backflow for values of the a 
less than some critical value. This result and feature 
was confirmed by Regev & Gitelman (2002), who 
abandoned the polytropic assumption and included 
an energy equation (employing the diffusion approx- 
imation in the treatment of the vertical radiative 
transport), showing that the polytropic assumption 
makes only very little substantive difference from 
the steady meridional flow solution of KK. 

The next step of the procedure consists of expand- 
ing all dependent variables asymptotic series in e, 
for example 



(6) 



Similar expansions are used for p and v, while the 
one for u consists of only odd powers of e. The fea- 
sibility of this choice is justified in Umurhan et al. 
(2006) (hereafter UNRS). These expansions are then 
substituted into the full equations HE]) and the ex- 
pression in various orders in e are collected. The re- 
sulting equations are similar to the ones found by 
KK, but there are some important differences: 

(i) Time-dependence is included so as to be able 
to study the dynamics and consequently time- 
derivatives of the dependent variables appear 
(the time unit taken as the typical rotation 
time, t = 1/fi). The problem is thus formu- 
lated as an initial value problem. 

(ii) Time-dependence is excluded in the lowest or- 
der term in the expansions and it is introduced 
only in the next significant order, in the form 
of an additive term, denoted by a prime, and 
can be viewed as a perturbation on the steady- 
state, but not necessarily an infinitesimal one) 



Thus the typical even and odd order expan- 
sions are of the form 



p(r,z,t) = p a {r, z) + e 2 [p 2 (r, z) + p' 2 (r,z,t)] 
+ e 4 [p 4 (r,z)+^(r,z, <)] + ••• (7) 
u(r, z,t) = e [ui(r, z) + u[(r, z, t)] + 

+ e 3 [u 3 (r,z) + u' 3 (r,z,t)} + --- (8) 

(iii) As in KK the solution of O (1) equations yield 
a steady disk of the Shakura & Sunyaev kind, 
but for a polytrope (see also UNRS for the 
case of an arbitrary polytropic index n) . The 
O (e) 2 equations can be conveniently split 
into a steady and time-dependent parts, ow- 
ing to the additive construction of the time- 
dependent disturbances. The steady part 
yields the KK steady solution, with possible 
backflows (which occur for a < a c ~ 0.7). 
We conclude this section by formulating the ini- 
tial value problem resulting from the O (e) 2 time- 
dependent equations. 



dtp' 2 



-drirpou'^ - d z (p v' 2 ), 



1 



Po 

4 1 

d t v 2 = -d z m 2 + - — d z (r)d z v' 2 ) H d z (r)rd z u x ) 

3po Por 



(9) 
(10) 



3po 



-d r (ru{] 
r 



d t n' 2 = - 



r z dr po 



(11) 
(12) 



where the zero-indexed quantities are known from 
the steady O (1) (Shakura-Sunyaev disk) solution 
(e.g., Slo = r~ 3 / 2 ) and we base the expression for r\ 
on the 0(1) solution only, which using ^ is 77 = 
(2/3)a Pq +1 ^" r~ 3 / 2 . In addition, from the poly- 
tropic law it follows that 1372 = wop' 2 / pq. 

The initial value problem defined by (|9l fT2)) must 
be complemented by boundary conditions. We have 
opted for retaining the integral steady mass inflow 
condition as in KK. In addition, guided by the phys- 
ical consideration that there should be no work done 
on the disk surface, we have required that the ex- 
ternal pressure and viscous stresses vanish at this 
surface. For details see UNRS. 
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3. Representative solutions 

Consider now solutions to the initial value prob- 
lem, posed in the previous section and for the sake of 
simplicity we set, from here and on, the polytropic 
index to the value n = 3/2 (as in KK). All solutions 
obviously depend on the initial conditions and our 
purpose here is to discuss some relevant represen- 
tative solutions, that can exhibit TG. It is impor- 
tant to note that the above four equations support 
two types of solutions, namely a pair comprising of 
purely vertical disturbances and a second pair con- 
taining disturbances in all the velocity components. 
This follows from the fact that equations (TIT)]) and 
(fT2| dynamically decouple from the other two, i.e., 
the quantities and tt' 2 are not dynamically influ- 
enced by v 2 and p' 2 , while the converse is not true. 
Thus, in principle, if perturbations are chosen such 
that u\ = Q' 2 — initially, they will remain so at 
all times and the motion will be of purely vertical 
acoustics (VA), as embodied in v 2 and p' 2 . On the 
other hand, if one chooses non-vanishing & f2' 2 
as initial conditions, the resulting dynamics com- 
prises of all variables. We shall call these solutions 
"driven" general acoustics (GA). This mathemati- 
cal structure is reminiscent of the behavior of the 
steady-state equations of KK and we shall now elab- 
orate on these two solution types. 

3.1. Meridional velocity disturbances 

It is quite straightforward to derive from the dif- 
ferential set (f9l ll~2"Tl two second order equations for 
just the function u[(r, z,t) and v' 2 {r,z,t), that is, 
the lowest order time-dependent meridional compo- 
nents of the velocity perturbations. They are 



Vu\ = 0, 

£v 2 = [d t F + g] u[, 



(13) 
(14) 



where the differential operators are defined (as op- 
erating on an arbitrary function ip(r, z,t)) in the fol- 
lowing way 



V = 



t - —d z (nd z ) 
Po 



d t ip- —d z (nd z ip) 

Po 



2 5 
C = dfip - -ru d z ip - -d z (m d z ip) - (9fro ) if 



g = d z 



■^—dJridJrip)] + —d r (rr]d z ip), 
3po rp 



2wq 
3rp 



d r (rp tp) 



(15) 



The boundary conditions completing this system 
can be specified using the 'free boundary' prescrip- 
tion at the disk surface, as mentioned at the end of 
the previous section (see UNRS for an explicit for- 
mulation). After solutions to equations (|13H14[) are 
found, it is possible to return to the original set of 
equations in order to find the two additional un- 
known functions p' 2 {r, z, t) and Q' 2 (r, z,t). 

We were able to find analytical solutions to the 
system (|13H14[) . In describing them we shall first dis- 
cuss, in subsection 13.21 the case u[ = and fl' 2 =0, 
which result from the special initial conditions in 
which these perturbations are not excited. This will 
thus allow us to be limited to solving only the homo- 
geneous part of (|14jl . The solutions of this equation 
constitute the above mentioned first pair of modes - 
VA. The solutions to the full set (|13H14j) in the gen- 
eral (i.e. inhomogeneous case, with u[ ^ 0) will be 
discussed in subsection 13.31 This pair of modes was 
called before GA and we remark once again that 
their dynamics can be viewed as being driven (by 
the inhomogeneous part). Finally, in subsection !3.4l 
we shall describe in some detail the TG behavior 
found in our solutions. A full exposition of the calcu- 
lations leading to the solutions quoted and discussed 
in subsections 13 . 2ll3~4l is presented in UNRS and its 
Appendices. In the following three subsections only 
the salient features will be given. 

3.2. Vertical acoustics (VA) 

The general solution to the inhomogeneous equa- 
tion p4p can be written, in the form 



v'i = v' h + v' 



(16) 



3po 



d z [r]d z (d t <p)}, 



where the indexes h and p stand for a solution of 
the homogenous equation and a particular solution 
of the inhomogeneous equation, respectively. As we 
have already remarked, the thinness of the disk al- 
lows for the density and vertical velocity fluctuations 
not to induce either radial or azimuthal motions at 
these orders. This implies that homogeneous solu- 
tions of ([H]), (i.e. with the equation's RHS set to 
zero, because — Q' 2 = say, see above) are per- 
fectly acceptable. 

These VA solutions have the form (see UNRS) 



5 



v 'h = v h (z, r) exp(cr v T) + c.c. with T = -^,(17) 

where the spatial eigenfunctions v' h are composed 
of the associated Legendre functions. This form of 
v' h already indicates that these solutions are insep- 
arable in r and t. The eigenvalues, which appear in 
complex conjugate pairs, are functions of a and in 
general they also depend, and quite sensitively, on 
the index n (see UNRS), but as mentioned before 
we concentrate here on the case n = 3/2 only. The 
eigenfrequency of the fundamental mode, for our 
chosen polytropic index is given by 



o\, = a ± i 

9 



1G 
81' 



1/2 



(18) 



All such eigenfrequencies have a negative real part 
and thus, as can be expected, the fundamental as 
well as all overtones, show temporal decay. 

3.3. General acoustics (G A) 

We turn now to the general solution to (|13M14p 
by first focusing on (|13|) . This equation admits an 
infinite set of cigenmode solutions in a similar way 
to the homogeneous solutions discussed before. The 
general derivation and structure of these eigenmodes 
is detailed in UNRS. For the sake of clarity, the dis- 
cussion here focuses only upon the dynamics associ- 
ated with the fundamental mode of the operator V '. 
This solution is given by (see UNRS) 



i'i = u'(z, r) exp(a u T + i?) + c.c. 



(19) 



with T = i/r 3 / 2 , as before, and in which the spa- 
tial eigenfunction u'{z, r) has the particularly sim- 
ple structure 



u'{r, z) = A(r) 



z 

1? 



(20) 



where h, the disk thickness is a known function of 
r, see (|22|) . A(r) is an amplitude whose radial func- 
tional form is technically arbitrary and would be set 
by the initial condition. An arbitrary phase factor 
of -d has also been introduced. The eigenvalues, de- 
noted in this case by o~ u , again come in complex con- 
jugate pairs and the fundamental frequency is given 

by 



la ± i. 



(21) 



Without loss of generality we may only consider the 
'+' solution since the phase factor $ can account 



for the '— ' solution. The temporal behavior of these 
modes is again one of decaying oscillations, but with 
frequencies given by the local rotation rate of the 
disk. Note that in the prescription for h(r) as given 
in UNRS for values of r sufficiently greater than the 
zero-torque radius, r„, (see KK) the disk height with 
respect to r is well approximated by, 



h(r) ~ hir, and thus 



dh 
dr 



hi, 



(22) 



where hi is determined by the mass flux rate and 
the value of a. Since hi simply adds a multiplicative 
factor to all the dynamical quantities, it plays no role 
in the quality of the ensuing evolution and, as such, 
we set hi = 1 without any loss of general flavor. 

With this solution for u[ we may find a particular 
solution to the vertical velocity, v' p by solving (fTi]) , 
together with (J9j) directly. The details of this, rather 
lengthy procedure are presented in UNRS, but we 
can highlight the major steps here. We posit the 
following Ansatz for v' p and p' p , 

v' p = v p ((,r)e^ T + V p (( 1 r)(Te°» T ) + c.c, (23) 
p' p = p p ((,r)e°» T + R p ((,r)(Te°» T ) + c.c, (24) 

where £ = z/h(r). 

Recall that we solve only for the fundamental 
mode of the system (see UNRS) and that there is ob- 
viously a dependence on a. It is important to notice 
the explicit appearance of a multiplicative T term 
in the above expressions. This temporal functional 
dependence is necessary in order to balance terms 
arising on the RHS of (fT4|) - the 'driving' term and 
will, of course, cause TG (see below in the next sub- 
section). In general, the form of these lowest order 
structure eigenfunctions is 



V„ 




(25) 



namely, it is a polynomial in £ of the third degree 
with only odd powers in £. The four coefficients are 
functions of r and a, i.e. a t = ai(r,a) and bi = 
bi(r, a). 

The solution to the density perturbation, p' p can 
be quite straightforwardly obtained using (|20|). (|23|) 
and , taking the functional form of po as detailed 
in the O (1) steady-state solution (see KK, UNRS) 
and finally integrating equation ([9]) with respect to 
time. This gives 



G 



Pp 



(w 2 ) 



1 

2\ 2 




,(26) 



where, again, the coefficients are functions of r and 
a, that is, Cj = Ci(r, a) and di = di(r, a). 

The polynomials appearing in the square brack- 
ets of the expressions for p p and R p have only even 
powers of £. The detailed forms of the coefficients 
depend on the form of the generalized initial pertur- 
bation radial structure A(r) (as well as on n, in the 
general polytrope case). We avoid presenting them 
here because they are very long and cumbersome; 
and the details of the coefficients have no effect on 
the resulting transient dynamics, since they describe 
only a particular type of disturbance. For the sake 
of simplicity most analytic results which will be pre- 
sented below assume A = e J7r / 4 (that is, a partic- 
ularly simple form of the initial conditions with no 
dependence of the disturbance). We shall however 
show also one instance of a more complicated initial 
condition consisting of a localized Gaussian profile 



3.4. Transient dynamics and growth 



As is evident by inspection of ([23]) . v' p and p' p ex- 
hibit a pronounced TG during a finite time interval 
(due to the factor T multiplying the decaying ex- 
ponential), before eventually decaying. To demon- 
strate this more clearly we may consider some in- 
tegrated energy quantities. The first of these is an 
energy density, per unit surface of the disk 



£ (r,i; a) 




(27) 



£ a (r, t; a) is to be interpreted as the acoustic energy 
(per unit area of the disk) consisting of the kinetic 
energy in the vertical velocity disturbances and the 
compression energy due to the density disturbances. 
The contribution of the (purely oscillatory) velocity 
components, resulting from the homogeneous part 
of the disturbance equations, are left out of (|27|) . 

£ a depends on the a parameter as well as the 
structure of the radial velocity perturbation A{r) 
(cf. [20|) . The choice of A(r) as constant with respect 
to r makes the integral defined in (|27|) analytically 



Comparative TP: Differing a at r = 1 




T = tr 



Fig. 1. The TG in the quantity £ a for the fundamental mode 
at r = 1, n = 3/2 and where A{r) = e i7r / 4 , dA/dr = 0. 
The four curves correspond to a = 0.1,0.01,0.001 and are 
presented on a log-log plot. All curves are £ a , scaled to its 
value at T = 0. The rise time, as well as the maximal value 
achieved, are proportional to the inverse of a. £ a also exhibits 
oscillations with period it that sit atop the general transient 
trend. 

tractable. Nonetheless the expressions, which have 
been verified with the aid of Mathematica 5.0, are 
very long and will not be displayed here; only their 
essential features will be addressed. With A = e l7r / 4 , 



£ a has the functional form 



,-2aT 



£ a (r,t;a) 



r 3/2 



-T(T, cos 2T, sin 2T; a), 



(28) 



where T is a well-defined analytical, albeit compli- 
cated, function of its arguments. Thus, £ a depends 
on time only through the variable T — t/r 3 / 2 , which 
is the time measured in units of the local disk rota- 
tion period, at r, divided by 2-7T. 

In Figure [TJ the evolution of £ a as a function of 
the variable T is shown for different values of a. T 
is actually a similarity variable and the value of r 
explicitly appears only in the pre-factor multiplying 
the function T. For the case displayed in the figure 
the radius is fixed (set to r = 1). Different values 
of r will merely change the overall height of the re- 
sponse (oc r~( 3 / 2 )) while the shape of the function 
is self-similar and is always decaying at long times 
T. Inspection of the figure, uncovers a very impor- 
tant property of the TG. As a decreases, the magni- 
tude of the growth becomes more promineniPl (for 
example, it is up to ~ 3 orders of magnitude for a = 
0.001) and the maximum occurs at correspondingly 



1 For an inviscid disk, i.e. with formally a = 0, Umurhan 
& Shaviv (2006) obtained algebraic growth with no decay, 
which is quite telling, but should be considered carefully (see 
Sec. [J] on the validity of asymptotic expansions). 
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Disturbance Energy 
a = 0.01 , e = 0.05 








□ 
D 
P 
□ 



r: radius (AU) 



Fig. 2. TG in the vertical velocity disturbance kinetic energy 
density in a disk with a = 0.01, e = 0.05. A(r) is as in l|30|l 
with ro = 3 and A = 2.5. Time is given in terms of disk 
rotation times at the radius r = 1, so that T = t. The mag- 
nitude of the disturbance energy density is indicated by the 
shadings (whose meaning is indicated on the right of each 
panel) on a cut through the disk's meridional plane. Signif- 
icant growth in the disturbance energy is observed (almost 
three orders of magnitude). Crenellation patterns that ap- 
pear in the fluctuating boundary of the disk are ultimately 
wiped out by viscous decay. 

later times. The time corresponding to the maxi- 
mum amplitude T max is roughly ~ 5/ (8a) in thij^l 
case. The rise in energy is modulated by oscillations 
which arise from the fact that there are correlations 
between pairs of variables which contribute to the 
overall growth in £ a (T). These variables oscillate in 
a frequency defined by the imaginary part of the 
eigenvalue (f2~Tj) and because £ a (T) is composed of 
products of pairs of dynamical variables terms de- 
pending on 2T appear - explaining why the observed 
period is half the orbital period at a given radius. 

£ a is a function of the radial position r, but one 
may form the integral quantity E a (t;a), i.e., the 
total disturbance acoustic energy contained in a 
global portion of the disk (a ring), by integrating 
over the radial range in which these disturbances 
are assumed to exist, that is, between the inner and 
outer bounds r m i n and r max of the ring. 



E a (t;a) 



£ a (r, t] a) 2nr dr. 



(29) 



2 For the fe-th overtone, T m 
1/k 2 , see UNRS) 



will be smaller by a factor of 



We have found (see UNRS) that A(r) affects the 
spatial details but not the global time behavior and 
that the occurrence of the peak energy value is de- 
layed for larger values of r max , at given values of a 
(and r m i„), although the value of the energy at the 
peak time remains roughly the same for a given a. 
To get an idea about the expected qualitative spatio- 
temporal behavior of such transiently growing so- 
lutions we have computed the results for a case in 
which A(r) is not constant (and, thus, the possibil- 
ity to evaluate the relevant integrals analytically can 
not be expected, in general). We consider here an 
initial condition in which the radial velocity ampli- 
tude has a Gaussian form around some fixed radius, 
r say, and having a width A, such that it is well 



contained within the region r m j n < r < r max 



A(r) 



e l7r/4 exp 



(r - r ) 2 



A 



/ 



(30) 



The fluid is assumed to be otherwise initially undis- 
turbed. As pointed out above, with this form for A 
we lose the ability to find an analytic solution and 
must resort to numerical evaluation. In Figure[5]we 
show, in a contour plot in the disk meridional plane, 
the kinetic energy (per unit volume) contained in the 
vertical velocity disturbance v p . The TG and decay 
of the disturbance is shown in the time sequence of 
figures displaying the spatial structure of the distur- 
bance. We also depict how the disk surface moves in 
response to the imposed perturbation by solving the 
equation of motion for the boundary at this order. 
Time is given in units of rotation times of the disk, 
as measured at the radius r = 1. Because the re- 
sponse of the disk surface is an integral with respect 
to time, it is not a surprise to see a large amplitude 
in the surface position long after the kinetic energy 
has started to die away. 

4. Summary and Discussion 

We conclude with some remarks on possible im- 
provements to the asymptotic analysis of the sort 
done here and prospects for the future, e.g., exten- 
sions to non-axisymmetric perturbations. Asymp- 
totic expansions, when viable, are often very robust 
and provide a good approximation to the solution 
when truncation to only few leading terms is done. 
Obviously, when a term in the series becomes very 
large it may 'break its order', that is, become larger 
than a previous term and as such make the expan- 
sion invalid in this region. In our expansions succes- 
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sive terms ratios are of O (e 2 ) , and thus the proce- 
dure's validity should not severely be limited even 
up to a growth factor of ~ 1000 or so (in the velocity 
or density perturbations) . This matter is further dis- 
cussed in UNRS, here we just remark that the valid- 
ity with respect to time for weakly viscous solutions 
are somewhat influenced by the thinness of the disk: 
smaller values of e mean that the solutions are valid 
for longer times after the initial disturbance. More 
importantly, however, for a given value of e one must 
not be too zealous or overreaching by attempting to 
infer the quantitative behavior of the disk for arbi- 
trarily small values of a - which, as one will recall, 
is formally assumed here to be an O (1) quantity. 

Despite these caveats, the procedure when car- 
ried to higher order introduces corrections which are 
technically non-linear. Careful consideration must 
be undertaken in order to handle the response at 
these higher orders. This may entail treating the dis- 
turbance amplitudes for the lower order solutions 
(like A(r) in u[) as weakly non-linear governed by 
a second 'slow' time (e.g. the amplitude is instead 
written as A = A(r,r) where r = e 2 t) in a man- 
ner analogous to the treatment of non-linear thick 
polytropes (e.g, Balmforth & Spiegel, 1996). 

The approach used here may be generalized 
in a number of directions. Allowing for non- 
axisymmctric perturbations, including the disk in- 
ner and outer boundary in some kind of boundary 
layer analysis and relaxing the polytropic assump- 
tion seem to be the most obvious generalizations. 
We have found the presence of prominent TG in 
the simplest cases. It is difficult to imagine that it 
will be suppressed in the more general conditions 
although the effect of radiative energy losses on TG 
must be carefully examined. 

The question concerning the ultimate fate of the 
transiently grown perturbations and their ability to 
induce a state of sustained complex dynamical ac- 
tivity in the disk remains open. In this context it is 
worthwhile to notice that since the TG decay times 
are of the order of hundreds of rotation periods, it is 
conceivable that AD, which are usually not isolated 
systems, may experience recurrent external pertur- 
bations on such time scales and in this way the dy- 
namical activity may be sustained. Extensive nu- 
merical calculations of AD are however needed to 
decide if TG may lead, through non-linear processes, 
to sustained turbulence, or at least a dynamical 
state in which adequate angular momentum trans- 
port can be sustained. Such high-resolution global 
3-D calculations are, however, still above the abil- 



ity of the present computer power and it may be 
thus advantageous to also exploit sophisticated non- 
linear asymptotic methods to complement and guide 
them. 

5. Afterword 

Accretion disks were introduced to explain astro- 
nomical X-ray sources which were discovered close to 
four decades ago, when Jean-Pierre Lasota was still 
a graduate student. At the time of the conference 
held at the Trzebieszowice Castle, on the occasion 
of his 65-th birthday, it seemed quite obvious that 
although we have gone a long way towards under- 
standing the physics of these fascinating and astro- 
nomically ubiquitous objects, much still remains to 
be elucidated. In order to enable accretion rates that 
are compatible with observations of such diverse ob- 
jects as forming stars, close binary systems or ac- 
tive galactic nuclei, an efficient physical mechanism 
for angular momentum transport must be present in 
these AD flows. It has been recognized at the outset 
that fluid turbulence offers a natural mechanism of 
this sort, but hydrodynamical (and also MHD) tur- 
bulence remains, unfortunately, until this day to be 
one of the major bugbears of theoretical physics. 

Shear flows are known to be particularly difficult 
in trying to understand the physics of transition to 
turbulence in them. Even the simplest flows of this 
kind, like the plane Couette and Poiseulle ones and 
the pipe flow, are not well understood in this respect 
(see, e.g., Drazin 2002). How, then, can we reason- 
ably expect to adequately understand this problem 
in AD, with their literally astronomical Re, strong 
shear and rotation, compressibility, magnetic fields 
permeating an ionized medium and other physical 
complications? The lesson learned from stellar tur- 
bulent convection is that the crudest (e.g, mixing- 
length) phcnomenological approaches may be very 
effective in constructing observationally viable mod- 
els of stars, but detailed understanding is extremely 
difficult. 

As far as AD are concerned, making significant 
progress beyond the famous a model requires, at 
least in my view, a combined synergetic effort of 
different groups of researchers, approaching the 
problem from many different angles and employ- 
ing essentially all relevant methods that theoretical 
physics can offer. Phenomenology, analytical and 
semi-analytical investigation of simplified systems, 
laboratory experiments and numerical simulations 
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of various kinds can all be useful, and especially 
if they fertilize each other. Paradigms have always 
played an important role in scientific progress, of- 
ten consecutively replacing one another. Dogmas, 
however, have always been stumbling blocks, es- 
pecially in the education and mentoring of young 
researches who are naturally able to offer fresh and 
independent ideas. 

If we are able to avoid such dogmas and truly co- 
operate in the way hinted on above, there are good 
chances that the physics of accretion disks will be 
much better understood when we celebrate Jean- 
Pierre's proverbial 120-th birthday, and hopefully 
even much sooner. 
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